Benchmarking Cross-Docking Strategies for Structure-Informed Machine Learning in Kinase Drug Discovery.

In recent years machine learning has transformed many aspects of the drug discovery process including small molecule design for which the prediction of the bioactivity is an integral part. Leveraging structural information about the interactions between a small molecule and its protein target has great potential for downstream machine learning scoring approaches, but is fundamentally limited by the accuracy with which protein:ligand complex structures can be predicted in a reliable and automated fashion. With the goal of finding practical approaches to generating useful kinase:inhibitor complex geometries for downstream machine learning scoring approaches, we present a kinase-centric docking benchmark assessing the performance of different classes of docking and pose selection strategies to assess how well experimentally observed binding modes are recapitulated in a realistic cross-docking scenario. The assembled benchmark data set focuses on the well-studied protein kinase family and comprises a subset of 589 protein structures co-crystallized with 423 ATP-competitive ligands. We find that the docking methods biased by the co-crystallized ligand-utilizing shape overlap with or without maximum common substructure matching-are more successful in recovering binding poses than standard physics-based docking alone. Also, docking into multiple structures significantly increases the chance to generate a low RMSD docking pose. Docking utilizing an approach that combines all three methods (Posit) into structures with the most similar co-crystallized ligands according to shape and electrostatics proofed to be the most efficient way to reproduce binding poses achieving a success rate of 66.9 % across all included systems. The studied docking and pose selection strategies-which utilize the OpenEye Toolkit-were implemented into pipelines of the KinoML framework allowing automated and reliable protein:ligand complex generation for future downstream machine learning tasks. Although focused on protein kinases, we believe the general findings can also be transferred to other protein families.

Prioritizing Small Sets of Molecules for Synthesis through in-silico Tools: A Comparison of Common Ranking Methods.

Prioritizing molecules for synthesis is a key role of computational methods within medicinal chemistry. Multiple tools exist for ranking molecules, from the cheap and popular molecular docking methods to more computationally expensive molecular-dynamics (MD)-based methods. It is often questioned whether the accuracy of the more rigorous methods justifies the higher computational cost and associated calculation time. Here, we compared the performance on ranking the binding of small molecules for seven scoring functions from five docking programs, one end-point method (MM/GBSA), and two MD-based free energy methods (PMX, FEP+). We investigated 16 pharmaceutically relevant targets with a total of 423 known binders. The performance of docking methods for ligand ranking was strongly system dependent. We observed that MD-based methods predominantly outperformed docking algorithms and MM/GBSA calculations. Based on our results, we recommend the application of MD-based free energy methods for prioritization of molecules for synthesis in lead optimization, whenever feasible.

Discovery and Characterization of the Potent and Highly Selective 1,7-Naphthyridine-Based Inhibitors BAY-091 and BAY-297 of the Kinase PIP4K2A.

PIP4K2A is an insufficiently studied type II lipid kinase that catalyzes the conversion of phosphatidylinositol-5-phosphate (PI5P) into phosphatidylinositol 4,5-bisphosphate (PI4,5P2). The involvement of PIP4K2A/B in cancer has been suggested, particularly in the context of p53 mutant/null tumors. PIP4K2A/B depletion has been shown to induce tumor growth inhibition, possibly due to hyperactivation of AKT and reactive oxygen species-mediated apoptosis. Herein, we report the identification of the novel potent and highly selective inhibitors BAY-091 and BAY-297 of the kinase PIP4K2A by high-throughput screening and subsequent structure-based optimization. Cellular target engagement of BAY-091 and BAY-297 was demonstrated using cellular thermal shift assay technology. However, inhibition of PIP4K2A with BAY-091 or BAY-297 did not translate into the hypothesized mode of action and antiproliferative activity in p53-deficient tumor cells. Therefore, BAY-091 and BAY-297 serve as valuable chemical probes to study PIP4K2A signaling and its involvement in pathophysiological conditions such as cancer.

Discovery of the SMYD3 Inhibitor BAY-6035 Using Thermal Shift Assay (TSA)-Based High-Throughput Screening.

SMYD3 (SET and MYND domain-containing protein 3) is a protein lysine methyltransferase that was initially described as an H3K4 methyltransferase involved in transcriptional regulation. SMYD3 has been reported to methylate and regulate several nonhistone proteins relevant to cancer, including mitogen-activated protein kinase kinase kinase 2 (MAP3K2), vascular endothelial growth factor receptor 1 (VEGFR1), and the human epidermal growth factor receptor 2 (HER2). In addition, overexpression of SMYD3 has been linked to poor prognosis in certain cancers, suggesting SMYD3 as a potential oncogene and attractive cancer drug target. Here we report the discovery of a novel SMYD3 inhibitor. We performed a thermal shift assay (TSA)-based high-throughput screening (HTS) with 410,000 compounds and identified a novel benzodiazepine-based SMYD3 inhibitor series. Crystal structures revealed that this series binds to the substrate binding site and occupies the hydrophobic lysine binding pocket via an unprecedented hydrogen bonding pattern. Biochemical assays showed substrate competitive behavior. Following optimization and extensive biophysical validation with surface plasmon resonance (SPR) analysis and isothermal titration calorimetry (ITC), we identified BAY-6035, which shows nanomolar potency and selectivity against kinases and other PKMTs. Furthermore, BAY-6035 specifically inhibits methylation of MAP3K2 by SMYD3 in a cellular mechanistic assay with an IC50 <100 nM. Moreover, we describe a congeneric negative control to BAY-6035. In summary, BAY-6035 is a novel selective and potent SMYD3 inhibitor probe that will foster the exploration of the biological role of SMYD3 in diseased and nondiseased tissues.

The potent AMPK inhibitor BAY-3827 shows strong efficacy in androgen-dependent prostate cancer models.

PURPOSE: 5' adenosine monophosphate-activated kinase (AMPK) is an essential regulator of cellular energy homeostasis and has been associated with different pathologies, including cancer. Precisely defining the biological role of AMPK necessitates the availability of a potent and selective inhibitor. METHODS: High-throughput screening and chemical optimization were performed to identify a novel AMPK inhibitor. Cell proliferation and mechanistic assays, as well as gene expression analysis and chromatin immunoprecipitation were used to investigate the cellular impact as well as the crosstalk between lipid metabolism and androgen signaling in prostate cancer models. Also, fatty acid turnover was determined by examining lipid droplet formation. RESULTS: We identified BAY-3827 as a novel and potent AMPK inhibitor with additional activity against ribosomal 6 kinase (RSK) family members. It displays strong anti-proliferative effects in androgen-dependent prostate cancer cell lines. Analysis of genes involved in AMPK signaling revealed that the expression of those encoding 3-hydroxy-3-methyl-glutaryl-coenzyme A reductase (HMGCR), fatty acid synthase (FASN) and 6-phosphofructo-2-kinase/fructose-2,6-biphosphatase 2 (PFKFB2), all of which are involved in lipid metabolism, was strongly upregulated by androgen in responsive models. Chromatin immunoprecipitation DNA-sequencing (ChIP-seq) analysis identified several androgen receptor (AR) binding peaks in the HMGCR and PFKFB2 genes. BAY-3827 strongly down-regulated the expression of lipase E (LIPE), cAMP-dependent protein kinase type II-beta regulatory subunit (PRKAR2B) and serine-threonine kinase AKT3 in responsive prostate cancer cell lines. Also, the expression of members of the carnitine palmitoyl-transferase 1 (CPT1) family was inhibited by BAY-3827, and this was paralleled by impaired lipid flux. CONCLUSIONS: The availability of the potent inhibitor BAY-3827 will contribute to a better understanding of the role of AMPK signaling in cancer, especially in prostate cancer.

Discovery of a hidden transient state in all bromodomain families.

Bromodomains (BDs) are small protein modules that interact with acetylated marks in histones. These posttranslational modifications are pivotal to regulate gene expression, making BDs promising targets to treat several diseases. While the general structure of BDs is well known, their dynamical features and their interplay with other macromolecules are poorly understood, hampering the rational design of potent and selective inhibitors. Here, we combine extensive molecular dynamics simulations, Markov state modeling, and available structural data to reveal a transiently formed state that is conserved across all BD families. It involves the breaking of two backbone hydrogen bonds that anchor the ZA-loop with the αA helix, opening a cryptic pocket that partially occludes the one associated to histone binding. By analyzing more than 1,900 experimental structures, we unveil just two adopting the hidden state, explaining why it has been previously unnoticed and providing direct structural evidence for its existence. Our results suggest that this state is an allosteric regulatory switch for BDs, potentially related to a recently unveiled BD-DNA-binding mode.

Computationally Empowered Workflow Identifies Novel Covalent Allosteric Binders for KRASG12C.

Due to its frequent mutations in multiple lethal cancers, KRAS is one of the most-studied anticancer targets nowadays. Since the discovery of the druggable allosteric binding site containing a G12C mutation, KRASG12C has been the focus of attention in oncology research. We report here a computationally driven approach aimed at identifying novel and selective KRASG12C covalent inhibitors. The workflow involved initial enumeration of virtual molecules tailored for the KRAS allosteric binding site. Tools such as pharmacophore modeling, docking, and free-energy perturbations were deployed to prioritize the compounds with the best profiles. The synthesized naphthyridinone scaffold showed the ability to react with G12C and inhibit KRASG12C . Analogues were prepared to establish structure-activity relationships, while molecular dynamics simulations and crystallization of the inhibitor-KRASG12C complex highlighted an unprecedented binding mode.

Computational Macrocyclization: From de†novo Macrocycle Generation to Binding Affinity Estimation.

Macrocycles play an increasing role in drug discovery, but their synthesis is often demanding. Computational tools that suggest macrocyclization based on a known binding mode and that estimate the binding affinity of these macrocycles could have a substantial impact on the medicinal chemistry design process. For both tasks, we established a workflow with high practical value. For five diverse pharmaceutical targets we show that the effect of macrocyclization on binding can be calculated robustly and accurately. Applying this method to macrocycles designed by LigMac, a search tool for de novo macrocyclization, our results suggest that we have a robust protocol in hand to design macrocycles and prioritize them prior to synthesis.

Benzoisoquinolinediones as Potent and Selective Inhibitors of BRPF2 and TAF1/TAF1L Bromodomains.

Bromodomains (BD) are readers of lysine acetylation marks present in numerous proteins associated with chromatin. Here we describe a dual inhibitor of the bromodomain and PHD finger (BRPF) family member BRPF2 and the TATA box binding protein-associated factors TAF1 and TAF1L. These proteins are found in large chromatin complexes and play important roles in transcription regulation. The substituted benzoisoquinolinedione series was identified by high-throughput screening, and subsequent structure-activity relationship optimization allowed generation of low nanomolar BRPF2 BD inhibitors with strong selectivity against BRPF1 and BRPF3 BDs. In addition, a strong inhibition of TAF1/TAF1L BD2 was measured for most derivatives. The best compound of the series was BAY-299, which is a very potent, dual inhibitor with an IC50 of 67 nM for BRPF2 BD, 8 nM for TAF1 BD2, and 106 nM for TAF1L BD2. Importantly, no activity was measured for BRD4 BDs. Furthermore, cellular activity was evidenced using a BRPF2- or TAF1-histone H3.3 or H4 interaction assay.

Collaborating to improve the use of free-energy and other quantitative methods in drug discovery.

In May and August, 2016, several pharmaceutical companies convened to discuss and compare experiences with Free Energy Perturbation (FEP). This unusual synchronization of interest was prompted by Schrödinger's FEP+ implementation and offered the opportunity to share fresh studies with FEP and enable broader discussions on the topic. This article summarizes key conclusions of the meetings, including a path forward of actions for this group to aid the accelerated evaluation, application and development of free energy and related quantitative, structure-based design methods.

Discovery and Characterization of a Highly Potent and Selective Aminopyrazoline-Based in Vivo Probe (BAY-598) for the Protein Lysine Methyltransferase SMYD2.

Protein lysine methyltransferases have recently emerged as a new target class for the development of inhibitors that modulate gene transcription or signaling pathways. SET and MYND domain containing protein 2 (SMYD2) is a catalytic SET domain containing methyltransferase reported to monomethylate lysine residues on histone and nonhistone proteins. Although several studies have uncovered an important role of SMYD2 in promoting cancer by protein methylation, the biology of SMYD2 is far from being fully understood. Utilization of highly potent and selective chemical probes for target validation has emerged as a concept which circumvents possible limitations of knockdown experiments and, in particular, could result in an improved exploration of drug targets with a complex underlying biology. Here, we report the development of a potent, selective, and cell-active, substrate-competitive inhibitor of SMYD2, which is the first reported inhibitor suitable for in vivo target validation studies in rodents.

Accuracy assessment and automation of free energy calculations for drug design.

As the free energy of binding of a ligand to its target is one of the crucial optimization parameters in drug design, its accurate prediction is highly desirable. In the present study we have assessed the average accuracy of free energy calculations for a total of 92 ligands binding to five different targets. To make this study and future larger scale applications possible we automated the setup procedure. Starting from user defined binding modes, the procedure decides which ligands to connect via a perturbation based on maximum common substructure criteria and produces all necessary parameter files for free energy calculations in AMBER 11. For the systems investigated, errors due to insufficient sampling were found to be substantial in some cases whereas differences in estimators (thermodynamic integration (TI) versus multistate Bennett acceptance ratio (MBAR)) were found to be negligible. Analytical uncertainty estimates calculated from a single free energy calculation were found to be much smaller than the sample standard deviation obtained from two independent free energy calculations. Agreement with experiment was found to be system dependent ranging from excellent to mediocre (RMSE = [0.9, 8.2, 4.7, 5.7, 8.7] kJ/mol). When restricting analyses to free energy calculations with sample standard deviations below 1 kJ/mol agreement with experiment improved (RMSE = [0.8, 6.9, 1.8, 3.9, 5.6] kJ/mol).

Mining electronic laboratory notebooks: analysis, retrosynthesis, and reaction based enumeration.

An approach to automatically analyze and use the knowledge contained in electronic laboratory notebooks (ELNs) has been developed. Reactions were reduced to their reactive center and converted to a string representation (SMIRKS) which formed the basis for reaction classification and in silico (retro-)synthesis. Of the SMIRKS that occurred at least five times, 98% successfully regenerated the original product. The extracted reaction rules (SMIRKS) and corresponding reactants span a virtual chemical space which showed a strong dependence on the size of the reactive center. Whereas relatively few robust reaction types were sufficient to describe a large part of all reactions, considerably more reaction rules were necessary to cover all reactions in the ELN. Furthermore, reaction sequences were extracted to identify frequent combinations and diversifying reaction steps. Based on the extracted knowledge a (retro-)synthesis tool was built allowing for de novo design of compounds which have a high chance of being synthetically accessible. In an example application of the de novo design tool, various feasible retrosynthetic routes to the query molecule were obtained. Reaction based enumeration along the top ranked route yielded a library of 29 920 compounds with diverse properties, 99.9% of which are novel in the sense that they are unknown to the public domain.

Calculation of relative free energies for ligand-protein binding, solvation, and conformational transitions using the GROMOS software.

The calculation of the relative free energies of ligand-protein binding, of solvation for different compounds, and of different conformational states of a polypeptide is of considerable interest in the design or selection of potential enzyme inhibitors. Since such processes in aqueous solution generally comprise energetic and entropic contributions from many molecular configurations, adequate sampling of the relevant parts of configurational space is required and can be achieved through molecular dynamics simulations. Various techniques to obtain converged ensemble averages and their implementation in the GROMOS software for biomolecular simulation are discussed, and examples of their application to biomolecules in aqueous solution are given.

Comparison of enveloping distribution sampling and thermodynamic integration to calculate binding free energies of phenylethanolamine N-methyltransferase inhibitors.

The relative binding free energy between two ligands to a specific protein can be obtained using various computational methods. The more accurate and also computationally more demanding techniques are the so-called free energy methods which use conformational sampling from molecular dynamics or Monte Carlo simulations to generate thermodynamic averages. Two such widely applied methods are the thermodynamic integration (TI) and the recently introduced enveloping distribution sampling (EDS) methods. In both cases relative binding free energies are obtained through the alchemical perturbations of one ligand into another in water and inside the binding pocket of the protein. TI requires many separate simulations and the specification of a pathway along which the system is perturbed from one ligand to another. Using the EDS approach, only a single automatically derived reference state enveloping both end states needs to be sampled. In addition, the choice of an optimal pathway in TI calculations is not trivial and a poor choice may lead to poor convergence along the pathway. Given this, EDS is expected to be a valuable and computationally efficient alternative to TI. In this study, the performances of these two methods are compared using the binding of ten tetrahydroisoquinoline derivatives to phenylethanolamine N-transferase as an example. The ligands involve a diverse set of functional groups leading to a wide range of free energy differences. In addition, two different schemes to determine automatically the EDS reference state parameters and two different topology approaches are compared.

Basic ingredients of free energy calculations: a review.

Methods to compute free energy differences between different states of a molecular system are reviewed with the aim of identifying their basic ingredients and their utility when applied in practice to biomolecular systems. A free energy calculation is comprised of three basic components: (i) a suitable model or Hamiltonian, (ii) a sampling protocol with which one can generate a representative ensemble of molecular configurations, and (iii) an estimator of the free energy difference itself. Alternative sampling protocols can be distinguished according to whether one or more states are to be sampled. In cases where only a single state is considered, six alternative techniques could be distinguished: (i) changing the dynamics, (ii) deforming the energy surface, (iii) extending the dimensionality, (iv) perturbing the forces, (v) reducing the number of degrees of freedom, and (vi) multi-copy approaches. In cases where multiple states are to be sampled, the three primary techniques are staging, importance sampling, and adiabatic decoupling. Estimators of the free energy can be classified as global methods that either count the number of times a given state is sampled or use energy differences. Or, they can be classified as local methods that either make use of the force or are based on transition probabilities. Finally, this overview of the available techniques and how they can be best used in a practical context is aimed at helping the reader choose the most appropriate combination of approaches for the biomolecular system, Hamiltonian and free energy difference of interest.

Effective-mode representation of non-Markovian dynamics: a hierarchical approximation of the spectral density. II. Application to environment-induced nonadiabatic dynamics.

The non-Markovian approach developed in the companion paper [Hughes et al., J. Chem. Phys. 131, 024109 (2009)], which employs a hierarchical series of approximate spectral densities, is extended to the treatment of nonadiabatic dynamics of coupled electronic states. We focus on a spin-boson-type Hamiltonian including a subset of system vibrational modes which are treated without any approximation, while a set of bath modes is transformed to a chain of effective modes and treated in a reduced-dimensional space. Only the first member of the chain is coupled to the electronic subsystem. The chain construction can be truncated at successive orders and is terminated by a Markovian closure acting on the end of the chain. From this Mori-type construction, a hierarchy of approximate spectral densities is obtained which approach the true bath spectral density with increasing accuracy. Applications are presented for the dynamics of a vibronic subsystem comprising a high-frequency mode and interacting with a low-frequency bath. The bath is shown to have a striking effect on the nonadiabatic dynamics, which can be rationalized in the effective-mode picture. A reduced two-dimensional subspace is constructed which accounts for the essential features of the nonadiabatic process induced by the effective environmental mode. Electronic coherence is found to be preserved on the shortest time scale determined by the effective mode, while decoherence sets in on a longer time scale. Numerical simulations are carried out using either an explicit wave function representation of the system and overall bath or else an explicit representation of the system and effective-mode part in conjunction with a Caldeira-Leggett master equation.

Effective-mode representation of non-Markovian dynamics: a hierarchical approximation of the spectral density. I. Application to single surface dynamics.

An approach to non-Markovian system-environment dynamics is described which is based on the construction of a hierarchy of coupled effective environmental modes that is terminated by coupling the final member of the hierarchy to a Markovian bath. For an arbitrary environment, which is linearly coupled to the subsystem, the discretized spectral density is replaced by a series of approximate spectral densities involving an increasing number of effective modes. This series of approximants, which are constructed analytically in this paper, guarantees the accurate representation of the overall system-plus-bath dynamics up to increasing times. The hierarchical structure is manifested in the approximate spectral densities in the form of the imaginary part of a continued fraction similar to Mori theory. The results are described for cases where the hierarchy is truncated at the first-, second-, and third-order level. It is demonstrated that the results generated from a reduced density matrix equation of motion and large dimensional system-plus-bath wavepacket calculations are in excellent agreement. For the reduced density matrix calculations, the system and hierarchy of effective modes are treated explicitly and the effects of the bath on the final member of the hierarchy are described by the Caldeira-Leggett equation and its generalization to zero temperature.

Comparison of three enveloping distribution sampling Hamiltonians for the estimation of multiple free energy differences from a single simulation.

We test the performance of three different reference state Hamiltonians for enveloping distribution sampling (EDS). EDS is an implementation of umbrella sampling which allows estimation of various free energy differences "on the fly" from a single simulation. This is achieved by construction of a reference state, which envelopes the regions of configuration space important to the various end states of interest. The proposed Hamiltonians differ in the way energy barriers separating these regions of configuration space are reduced. The test system consisted of 17 disubstituted benzenes in water and in complex with alpha-cyclodextrin. The calculated free energy differences correlate with thermodynamic integration results (R(2) > 0.99 for the ligands in water and R(2) > 0.98 for the ligands in complex with alpha-cyclodextrin). One of the reference state Hamiltonians outperformed the others in sampling the configuration space important to the various end states. In this reference state not all barriers between all pairs of states are reduced. Instead a minimum spanning tree of states is calculated, which connects states that are "closest" in configuration space.

Simple, Efficient, and Reliable Computation of Multiple Free Energy Differences from a Single Simulation: A Reference Hamiltonian Parameter Update Scheme for Enveloping Distribution Sampling (EDS).

We present an automatic adaptive scheme which allows fast optimization of the reference Hamiltonian parameters in enveloping distribution sampling (EDS). Six different variants of the update scheme have been tested on a condensed phase test system which included the recurrent deletion and creation of complete water molecules in water. All six schemes gave accurate free energy estimates with absolute errors of up to 1 kJ/mol for the worst scheme and up to 0.1 kJ/mol for the best scheme. Configurational sampling is focused on the regions where the end state energy difference distributions intersect, explaining the high accuracy and precision of the free energy estimates. The new update scheme makes the application of EDS to other systems, e.g. in ligand binding studies, easy as no reference state Hamiltonian parameters have to be chosen by the user. The only necessary input are the Hamiltonians of the various end states involved.